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Abstract. - A dynamic model for failures in biological organisms is proposed and studied 
both analytically and numerically. Each cell in the organism becomes dead under sufficiently 
strong stress, and is then allowed to be healed with some probability, ft is found that unlike 
the case of no healing, the organism in general does not completely break down even in the 
presence of noise. Revealed is the characteristic time evolution that the system tends to resist 
the stress longer than the system without healing, followed by sudden breakdown with some 
fraction of cells surviving. When the noise is weak, the critical stress beyond which the system 
breaks down increases rapidly as the healing parameter is raised from zero, indicative of the 
importance of healing in biological systems. 



Many systems in nature under steady external driving are led into failures. Examples 
are diverse ranging from disordered heterogeneous materials, earthquakes, and even to social 
processes. The fiber bundle models, describing well these systems, have received considerable 
attention [1]. Most studies of fiber bundles are based on the recursive breaking dynamics 
at discrete time step. One of the typical features in the models is that local stress arising 
from external driving tends to produce avalanches of microscopic failures, resulting in sta- 
tionary macroscopic breakdown of the system. The main issue here is thus the life time of 
a fiber bundle and the avalanche size distribution under various conditions [2-4]. There are 
shortcomings, however, in discrete dynamics: It cannot describe the time evolution in real 
continuous time, for which a remedy has been suggested recently [5]. 

Failures are also observed in biological systems, as a consequence of disease for instance; 
these are different from other failures in the sense that there may exist healing, in reparative 
or regenerative ways. Here the evolution to the stationary state appears more important than 
the stationary state itself. 

In this Letter we extend our previous work on the dynamic model for failures in fiber 
bundles [5] to biological systems consisting of cells. Incorporating the effects of healing into 
the model in a manner similar to that used for neural networks [6], we derive equations 
of motion for biological organisms in the form of delay-differential equations. From this 
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formulation, we obtain the evolution equation for the average number of living cells and find 
that there generally remains a finite fraction of living cells in the stationary state, manifesting 
the healing effects. This persists in the presence of noise, in sharp contrast to the case of no 
healing. The system is then explored numerically, which shows that the presence of healing 
assists the system to resist the stress longer before abrupt breakdown into the stationary 
state as well as increases the critical stress of external load that brings about the breakdown. 
The resulting behavior is reminiscent of the time course of degenerative disease progression 
such as diabetes, Alzheimer's disease, and possibly AIDS. We thus believe our approach may 
serve as a prototype model for disease progression at the cellular level, for which only a 
phcnomenological description is available so far [7]. 

We consider an organism consisting of N cells under external stress characterized by load 
F. Each cell has its own tolerance and endures the stress below the tolerance, thus remaining 
alive. The cell may become dead, however, if the tolerance is exceeded. We assign "spin" 
variables to these in such a way that Sj = —1 (+1) for the ith cell alive (dead). The state of 
the organism is then described by the configuration of all the cells, s = (si,s 2 , s N ). The 
total number N_ of living cells is related with the average spin s = N~ x Sj via 

^- = E^ = y(i-^ (1) 

i=i 

and we are interested in how A_ evolves in time as well as its stationary value. The total 
stress on the ith cell can then be written in the form 

TH = f + Y, V « L T i > (2) 

3 

where / = F/N is the stress directly due to the external load and Vij represents the stress 
transferred from the jth cell (in case that it is dead). The death of the ith cell with tolerance 
hi is determined according to: 

rji < hi =4> Si = -1 
rji > hi => Si = +1, 

which, in terms of the local field Ei = (rji — hi)(l — s)/2, can be simplified as Sj-Ej > 0. This 
determines the stationary configuration at which the system eventually arrives. 

For a more realistic description of the time evolution, we also take into consideration the 
uncertainty ("noise") present in real situations, which may arise from imperfections, random 
variations, and other environmental influences. We thus begin with the conditional probability 
that the ith cell is dead at time t+St, given that it is alive at time t. For sufficiently small St, 
we may write [6] 

p(8i= + l,t+5t\8i=-l,t;B',t-t d ) = ^[1+tanh/?^], (3) 

where s' = (s' l7 s' 2 , . . . , s' N ) represents the configuration of the system at time t— td and E[ — 
{Vi ~ hi){l — s')/2 is the local field at time t—td- Note the two time scales td and t r here: t d 
denotes the time delay during which the stress is redistributed among cells while the refractory 
period t r sets the relaxation time (or life time). The "temperature" T = /3 _1 measures the 
width of the tolerance region of the cells or the noise level: In the noiseless limit (T = 0) the 
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factor (1 + tanh/3i5 4 ')/2 in eq. reduces to the step function 0(E' i ), yielding the stationary- 
state condition. We also assign the non-zero conditional probability of the ith cell being 
repaired (regenerated) given that it is dead at time t, according to 

p( Si = - 1, t+6t\si= + 1, t; s', t-t d ) = (4) 

to 

where to is the time necessary for cell regeneration. Equations © and (@J can be combined 
to give a general expression for the conditional probability 

p(-Si, t+5t\si, t; s', t-td) = w i (s l ; s', t~t d )St (5) 

with the transition rate 



/ \ 1 
Wi{Si]S ,t-td) = — 



1\ f l\ 1-8 



2J ' \ a -2 



Si H — -tanh/3#< 



(6) 



where the healing parameter a = t r /to measures the relative time scale of relaxation and 
regeneration of a cell. Here we point out that td, t r , to, and accordingly a may be complicated 
functions of system properties like the average spin s and others, which may be incorporated 
into the model. In this work we restrict ourselves to the simplest case of these parameters 
being fixed. 

The behavior of the organism is then governed by the master equation, which describes 
the evolution of the joint probability P(s, t; s', t—tj) that the system is in state s' at time t—td 
and in state s at time t. Following the procedures in Ref. [5], we rescale time t in units of the 
delay time td and write accordingly the transition rate iOj(sj;s') = tdWi(si]s' ,t—td), which 
yields, in the limit 8t — > 0, 

-|p(M; s',i-l) = -J2[ w iisi;s')P{s,t;s',t-l)-w i (-s i ;s')P(F i s,t;s',t-l)} (7) 

i 

with FiS = (si, S2, — Sj, Sj+i, sjsr)- Here it has been noted that contributions from the 
intermediate configurations s' differing from s by only one cell survive; those from other 
configurations are of order (St) 2 or higher and thus vanish in the limit St — > 0. Then equations 
describing the time evolution of relevant physical quantities, with the average taken over 
P(s, t; s', t— 1), in general assume the form of differential-difference equations due to the delay 
in the stress redistribution. In particular, the average spin for the fcth cell, TOfe(t) = (s&) = 
J2 s ,s> s kP(s, t; s', t— 1) satisfies 

T Jt mk= Q~ a ) ~ Q+ a ) m fe+(^^ tanl1 ^) ( g ) 

where r = t r /t d gives the relaxation time (in units of td)- 

To proceed further, we assume equal load sharing, i.e., that the stress is distributed to 
every cell uniformly. In this case we have Vij — r\j /N- and accordingly, (1 — s)rji = 2/ from 
cq. . The infinite-range nature of equal load-sharing allows one to replace E' k by its average 
(E' k ) = f — (/ifc/2)[l — m(t— 1)], where it has been noted that s' is the configuration at time 
t—1, i.e., (s') = A^ -1 Ylj(s'j) = J2j m j{t~ 1) = rn(t—l). For convenience, we now rewrite 
eq. (jHJ in terms of the average number of living cells at time t. Defining Xk = (1 — rrik)/2 and 
x = N" 1 x k — (1 — ™)/2, we have, from eq. QJ, (N-) — Nx and thus obtain from eq. 
iJSJ the equation of motion for the average fraction of living cells: 

T j t x k(t) = ° - U + °) x * (*) - \ x &) tanh/3[/-/ijfex(t-l)]. (9) 
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We first examine the stationary solution of eq. with d,Xk/dt = 0: 

2a 

Xk ~ (l + 2a)+tanh/3(/-x/i fc )' ( ^ 

which, upon averaging over the tolerance distribution g(h), leads to the self-consistency equa- 
tion for the average fraction of living cells or the "health status" of the organism 

dhgih)- — ; r . (11) 

yy ' (l + 2a)+t&nhf3(f -xh) V ' 

It is obvious that regardless of noise, the complete breakdown, described by the null solution 
x — 0, is possible only for a — 0, i.e., when cells are not regenerated or repaired. This 
contrasts with the case a — 0, where noise usually gives rise to the breakdown of the system 
(x — 0) [5,8]. It is thus concluded that the biological systems become robust against noise 
due to the cell regeneration effects. 

In the noiseless limit (T = 0) we replace the term tanh/3(/— xh) in eq. (fTTf) by 9{f—xh) — 
8(xh—f ), to obtain 



dhg(h) 



o(hx-f) + — - e{f-hx) 

1 + a 



(12) 



For a continuous distribution of the tolerance, we formally solve eq. 11121) by simply performing 
the integration 

x = l- T ±-G(f/x), (13) 

where G{h) = jj 1 dhl g(h') is the cumulative distribution of the tolerance. Equation 1131) 
shows that a majority of cells in general becomes dead at stress above a critical value / c , 
which depends on a. This can be manifested with simple distributions, e.g., the bimodal 
distribution g(h) = pS(h-fi) + (1 — p)&(h—f2) with < f 2 < f±. A simple integration leads 
to 

1 for / < f 2 x 

{p + ^il + a)- 1 for f 2 x < f < fix 
a(l + a) -1 for fix < f. 

As / is raised from zero in the system with ji < pf\, the first breakdown occurs at / = 
yielding x = (p + a)(l + a) -1 , and subsequently the second one at / = (p + a)(l + a) _1 /i. 
For f 2 > pfi, on the other hand, there appears only one breakdown at / = /2. In these two 
cases the critical stress / c , beyond which the system breaks down significantly, is thus given by 
(p + a)(l + a)~ 1 fi and f 2 , respectively. Note that the survival fraction a(l + a) _1 is negligible 
for small healing parameter (a <C 1). 

We now turn our attention to the time evolution of the system. We consider first the 
simplest case of equal load sharing, and integrate directly the equation of motion in a 
numerical way. We consider the Gaussian distribution of the tolerance with unit mean and 
standard deviation a = 0.2, mostly in a system of N = 10 4 cells. Other distributions including 
the Weibull distribution have also been considered, only to give essentially the same results. 
Specifically, we have used ten different configurations of the tolerance distribution and set the 
relaxation time r = 50 and the time step At = 0.1. These parameter values have been varied, 
only to give no appreciable difference except for the time scale. 

We first display in fig. ^the health status of the organism at T = 0, with / being slightly 
larger than the critical stress f c {a=0) « 0.677. As reflected by the plateau, the system without 
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Fig. 1 - Evolution of the average fraction x of living cells with time t (in units of the delay time) 
for the Gaussian distribution of ft = 1 and a — 0.2 at T = 0. The stress is / = 0.685 slightly larger 
than / c (0) « 0.677. Solid, dotted, and broken lines correspond to a = 0, 0.03, and 0.06, respectively. 
Shown in the inset is the behavior for a — 0.6 (still at T — and a = 0). The stress / = 0.522 is 
again slightly larger than / c (0) (~ 0.518 for a = 0.6). 

healing resists the stress for some time before it breaks down completely. With small amount 
of healing, the duration of the plateau becomes longer and there remains a finite fraction of 
cells after the breakdown. We refer to this as the unhealthy state of the system [9]. Here 
further increase of a makes the duration of the plateau extremely long and the system does 
not break down, implying that the healthy state of the organism persists. For / is smaller 
than the critical stress / c (a=0), the average fraction x of living cells decreases very slowly (not 
shown) and the system remains healthy. Note also that the detailed degeneration behavior of 
the health status depends on the tolerance distribution. For example, a wider distribution in 
general brings about severer deterioration initially, followed by long duration of the unhealthy 
state before the breakdown (see the inset), similar to the progression of, e.g., AIDS. 

More interesting results are obtained at finite noise levels. For a = 0, non-negligible noise 
in the system eventually results in complete breakdown for any finite stress (i.e., f c = at 
T > 0.1 although there is numerical ambiguity at very low noise levels). When o / 0, on 
the other hand, the system may remain healthy up to some finite stress which depends on 
a, and thus f c does not vanish. To illustrate this, we plot in fig. HJa) the residual fraction 
x as a function of / for a = 0.1 at T = 0.1. It is clearly observed that the system becomes 
unhealthy abruptly near f c w 0.55 as / is increased. Figure H^b) displays how the critical 
stress f c varies with the healing parameter a at several noise levels. At T = 0, f c is observed 
to increase very slowly with a. In sharp contrast, at low but finite noise levels (T <C 1), which 
are presumably the case for biological systems in nature, f c increases very rapidly as a is 
raised from zero, revealing the crucial role of healing: The presence of even very weak healing 
can raise f c substantially from zero, thus assisting the system to resist moderate stress. At 
high noise levels, finite values of a are necessary to make f c nonzero. 

To further investigate this, we plot in fig. [3] the residual fraction Xf of living cells as a 
function of a for several values of / and T. Figure 01a) shows as a function of a for three 
different values of / near f c (a = 0) = 0.677 at T = 0. Observed is that sufficiently strong 
stress (i.e., large values of /) in general drives the system out of the healthy state. Otherwise, 
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Fig. 2 - (a) Residual fraction Xf of living cells versus / for a — 0.1 at T — 0.1. (b) Critical stress f c 
versus healing parameter a at several temperatures. 



as a is raised from zero for given value of /, xj grows from zero and increases sharply around 
the "critical value" a c depending on /, which is reminiscent of a phase transition, say, from 
the unhealthy to healthy state. The transition appears to become sharper in the presence of 
noise (T = 0.2) as shown in fig. |3Jb), with a c increasing with /. The transition, however, 




Fig. 3 - Residual fraction x / of living cells versus healing parameter a for several values of the initial 
stress / at noise level T = (a) 0, (b) 0.2, and (c) 0.4. Note the large error bars in the "transition 
region" in (a). Error bars in other regions are of the order of the symbol size. 
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tends to be continuous as the noise level is increased further [see fig. EI C )]- 

To confirm these findings, we have carried out Monte Carlo simulations directly and found 
that the results obtained from 20 independent runs with different initial configurations dis- 
play perfect agreement with those from numerical integrations of the same samples. In ad- 
dition, we have also examined the local load sharing case: Direct Monte Carlo simulations 
of two-dimensional systems, with the load being transferred to the shortest paths, show that 
qualitative behavior does not change except that the critical stress f c is smaller than that for 
the system with equal load sharing. Further, as the connectivity of each cell is increased, the 
behavior approaches that of the latter system [10]. This provides justification for the use of 
equal load sharing in a biological system, which may have enhanced connectivity in the two- 
or three-dimensional underlying structure. 

In summary, we have introduced a dynamic model for failures in biological organisms 
and investigated behaviors of the system under stress and healing. The dynamics takes into 
consideration uncertainty due to imperfections and environmental influences, described by 
noise. Such noise has been found to result in eventual breakdown of the system without 
healing. On the other hand, in the presence of healing, the system in general resists stress 
for a longer time compared with the case of no healing, and avoids complete breakdown. 
At weak noise the critical stress beyond which the system breaks down increases rapidly 
with the healing parameter, revealing the crucial role of healing. In spite of simplicity, the 
model shows quite interesting features, suggestive of the breakdown phenomena in biological 
systems. Finally, we point out that this model is rather general without details involved and 
thus provides a convenient starting point for wide potential applicability. In particular many 
refinements toward more realistic descriptions are conceivable. Corresponding applications to 
specific biological systems and their implications are left for further study. 
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